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Using black hole perturbation theory and arbitrary-precision computer algebra, we obtain the 
post-Newtonian (pN) expansions of the linear-in-mass-ratio corrections to the spin-precession angle 
and tidal invariants for a particle in circular orbit around a Schwarzschild black hole. We extract 
coefficients up to 20pN order from numerical results that are calculated with an accuracy greater than 
1 part in These results can be used to calibrate parameters in effective-one-body models of 

compact binaries, specifically the spin-orbit part of the effective Hamiltonian and the dynamically 
significant tidal part of the main radial potential of the effective metric. Our calculations are 
performed in a radiation gauge, which is known to be singular away from the particle. To overcome 
this irregularity, we define suitable Detweiler-Whiting singular and regular fields in this gauge, and 
we compute the invariants using mode-sum regularization in combination with averaging from two 
sides of the particle. The detailed justification of this regularization procedure will be presented in 
a forthcoming companion paper. 


I. INTRODUCTION 

Compact binaries are expected to be among the most 
significant sources for gravitational wave detectors, and 
over the last decade, modeling of these binaries has be¬ 
come a mature field. For two widely separated compact 
objects, post-Newtonian (pN) models of very high or¬ 
der have been developed mm- For compact objects of 
very dissimilar masses, calculations using gravitational 
self-force theory mm can now be performed with high 
accuracy at linear order in the mass ratio ISH3, and cal¬ 
culations at second order in the mass ratio are under¬ 
way [Mill- For objects of similar mass and small sepa¬ 
rations, fully nonlinear numerical relativity can be used, 
and it is continually pushing to cover more of the binary 
parameter space |15] . Furthermore, effective-one-body 
theory (FOB) [121 HZl has the potential to model all of 
these situations, given appropriate calibration data from 
each of the other models. 

In the last few years, work has increasingly been done 
at the interface between these various models [TH]. Of 
particular interest to us is that self-force data has been 
used to set benchmarks for numerical relativity, deter¬ 
mine high-order pN parameters, and calibrate FOB. 
One thing that has aided these efforts is the advent of 
extreme-accuracy numerical self-force results, with accu¬ 
racies of several hundred or thousands of digits mm- 
With these results, one can extract very high-order 
(linear-in-mass-ratio) post-Newtonian (pN) coefficients, 
which would be prohibitively difficult to calculate in an 
analytical pN framework. Perhaps surprisingly, with 
high-accuracy numerical self-force data one can extract 
not only numerical values of the pN coefficients, but their 
analytical form in terms of rational numbers (and their 
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product with tt), and for some terms, even analytical 
forms involving linear combinations of transcendentals. 
The details of this procedure are given in Ref. [20] . 

To use data from one model in another, one must calcu¬ 
late gauge- and coordinate-invariant quantities that can 
be translated between formalisms. These invariants can 
be divided into those related to dissipation and those re¬ 
lated to conservative effects. In the dissipative sector, for 
example, Ref. [TH] calculated the ingoing-at-horizon and 
outgoing-at-infinity gravitational-wave fluxes to 21pN or¬ 
der for a particle in circular, equatorial orbit around a 
Kerr black hole. However, most of the work done has 
been in the conservative sector, where a richer variety 
of invariants appear. The paragon of conservative invari¬ 
ants is Detweiler’s redshift invariant ED, which is related 
to the binding energy and angular momentum of a quasi¬ 
circular binary system. Reference [S] determined the pN 
expansion of this invariant up to I0.5pN order for a par¬ 
ticle in circular orbit around a Schwarzschild black hole, 
and Ref. [20] has pushed that expansion even further, to 
21.5pN order. 

Recently, other conservative invariants on circular or¬ 
bits in Schwarzschild spacetime have been discovered and 
calculated: the linear-in-mass-ratio correction to the spin 
precession (due to the “self-torque” effect) was numeri¬ 
cally calculated in Ref. [12] , where the pN expansion was 
found up to 3pN order (also see [13]); and conservative 
tidal invariants, the eigenvalues of certain electric and 
magnetic tidal tensors, were calculated in Ref. [8] to more 
than machine precision. 

In parallel with this numerical work, there have been 
analytical calculations of invariants, made possible by 
performing a pN expansion of the equations of self-force 
theory. Reassuringly, in all cases the results of the ana¬ 
lytical calculations have precisely agreed with the expan¬ 
sions extracted from numerical data. Reference [24] ana¬ 
lytically calculated the pN series of ingoing and outgoing 
fluxes for a particle in circular, equatorial orbit around a 
Kerr black hole and conhrmed the results in Ref. [H] to 
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llpN order. Bini & Damour have analytically calculated 
the pN expansion of the redshift invariant up to 9.5pN or¬ 
der [5SH5H], the spin precession to 10.5pN order [^150] . 
and tidal invariants to 7.5pN order m- 

Self-force calculations, whether numerical or analyt¬ 
ical, of these invariants have proven to be very effec¬ 
tive in calibrating EOB. Bini & Damour have used the 
spin-precession invariant to calculate the correspond¬ 
ing terms in the spin-orbit part of the effective Hamil¬ 
tonian through an effective gyro-gravitomagnetic ratio. 
They have likewise transcribed the tidal invariants into 
the tidal interaction energy in the EOB formalism. A 
semi-analytical simulation of tidally interacting binary 
neutron stars using these calibrations was performed in 
Ref. [32] , and for compact stars the model was confirmed 
by fully numerical simulation in Ref. |33j . 

Despite these successes, at this point in time, computa¬ 
tions of the spin-precession and tidal invariants have been 
performed at much lower accuracy than the extreme- 
accuracy calculations in Refs. PUS). In this paper, we 
pursue the program initiated in Refs. p in by calcu¬ 
lating the spin-precession and tidal invariants with an 
accuracy greater than 1 part in for circular orbits 

in Schwarzschild with orbital radii in the range 10^®- 
9 X This allows us to extract the pN expansion of 

the invariants to about 20 pN order. Such a procedure of 
extracting pN coefficients from high precision numerical 
data will be very important when calculating conserva¬ 
tive invariants for orbits in Kerr spacetime, where ana¬ 
lytical checks or computations are likely to be extremely 
difficult. 

Besides our high accuracy, our calculations differ from 
others in our choice of gauge. Dolan et al PI23 per¬ 
formed their calculations in the Lorenz gauge and Regge- 
Wheeler-Zerilli gauge. Bini & Damour I2S1I3I] performed 
theirs in the Regge-Wheeler-Zerilli gauge using analyt¬ 
ical solutions of the Regge-Wheeler equation given by 
Ref. [33]. In this paper we use a radiation gauge, fol¬ 
lowing the methods of Refs. PI35HS7I- Perturbations in 
radiation gauges can be reconstructed from solutions to 
the Teukolsky equation, making this gauge, unlike Lorenz 
and Regge-Wheeler-Zerilli, ideal for calculations in Kerr 
spacetime; hence, we view the calculations in this paper 
as a step toward computing analogous invariants in Kerr. 

Although it has this distinct advantage, the radiation 
gauge has historically had two drawbacks. First, the re¬ 
constructed metric does not contain the full solution of 
the linearized Einstein equation, but must be completed 
by adding a (non-radiation-gauge) perturbation that en¬ 
codes the particle’s mass and angular momentum. (We 
shall refer to the gauge of the completed solution as a 
completed radiation gauge; elsewhere, it has been called 
a modified radiation gauge.) Second, most quantities of 
interest in self-force theory are constructed from a cer¬ 
tain smooth vacuum perturbation called the Detweiler- 
Whiting (DW) regular field [38|, which in practice is 
obtained from the retarded metric perturbation by sub¬ 
tracting a certain divergent field called the DW singular 


field. Until recently, there has been no well-justified way 
of performing this procedure in a radiation gauge, which 
has long been known to have pathological singularities 
away from the particle [55] . 

However, both of these issues have, to a large extent, 
been dealt with. Regarding the completion problem, it 
has been incrementally better understood over the last 
few years [37||3nH131- The present paper provides some 
description of the requirements on the completion term 
in the metric perturbation in Schwarzschild, and a fu¬ 
ture paper will contain a fuller description of the prob¬ 
lem and its partial resolution in Kerr jUj . Regarding the 
problem of regularization, recent work of Pound, Merlin, 
and Barack (hereafter PMB) [33] has provided a rigor¬ 
ous means of obtaining the correct gravitational self-force 
from radiation-gauge perturbations, and this method has 
been concretely implemented in Ref. |45] . Here we use a 
regularization scheme, based on the results of PMB, that 
combines standard mode-sum regularization |46j with a 
certain averaging procedure. The rigorous justification 
of this scheme will be presented in a companion paper, 
but here we go some ways toward making it sensible by 
providing a suitable definition of the DW regular field in 
our completed radiation gauge. 

The article is organized as follows. In Sec. [TT] we de¬ 
scribe the setup of our calculations, the invariants we 
calculate, and our procedure for reconstructing and com¬ 
pleting the retarded metric perturbation. In Sec . [HI| we 
describe our regularization procedure. In Sec. |IV[ we 
discuss the extraction of pN coefficients. Most of these 
coefficients we determine as numerical values, but we also 
determine many coefficients in analytical form. Because 
of the length of the expansions, we present the analytical 
coefficients in an appendix]^ and we present the numer¬ 
ical coefficients in the accompanying supplementary files. 
We work in units with G = c = I. 

Before proceeding to the body of the paper, we note 
that in parallel with our efforts, a group comprising Chris 
Kavanagh, Adrian Ottewill and Barry Wardell [32] have 
also performed very high-order pN expansions of the 
same invariants, using analytical methods in the Regge- 
Wheeler-Zerilli gauge. At all orders that we have com¬ 
pared, these two independent calculations have yielded 
identical results. 


II. SETUP 

A. Circular orbits and invariant quantities 

We consider a particle of mass ^ orbiting a 
Schwarzschild black hole of mass M. At zeroth order 
in /r/M, the trajectory is a circular geodesic 70 of the 
Schwarzschild metric In Schwarzschild coordinates. 
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its 4-velocity is given by 


with 



( 1 ) 


where rg is the orbital radius, O is the orbital frequency, 
and and t°‘ are the rotational and timelike Killing 
vectors of the Schwarzschild spacetime, respectively. The 
orbital energy and angular momentum per unit mass are 
given by 


E = 


I _ 2 M 

■ro 



L, = 


/ 


Mrg 

^ _ 3M • 


( 2 ) 


At linear order in fi/M, the particle sources a pertur¬ 
bation that exerts a self-force on the particle. The 
perturbed trajectory 7 is no longer a geodesic of the back¬ 
ground metric, but it is a geodesic of an effective metric 
9in' = 9iiv + h^v, where is the DW regular field. If we 

consider only the conservative dynamics, the perturbed 
trajectory can be chosen, like the zeroth-order one, to be 
a circular orbit. 

The perturbed spacetime loses the background’s sta- 
tionarity and rotational symmetry, but it inherits the or¬ 
bit’s helical symmetry. In convenient gauges, the helical 
Killing vector is given by 


= r + 


(3) 


and at all points in the spacetime we have = 0 = 

^kh^u- We can use this Killing vector in combination 
with or derivatives of to construct meaningful 
scalar quantities on the orbit. Let be the four-velocity 
of the perturbed trajectory, normalized as = 

— 1. The simplest of the invariants is Detweiler’s redshift 
invariant. 


u* = (4) 

which is the rate of change of the time of an inertial frame 
at infinity relative to proper time in along the orbit. 
Next is the spin precession invariant, 

■0 = 1 — (5) 

which for a spinning particle is the angle of spin pre¬ 
cession per unit radian of the circular orbit. Here the 
precession frequency is where 

and is the covariant derivative com¬ 
patible with g^i^. The last set of quantities we shall 
consider are the eigenvalues of the electric and magnetic 
quadrupole tidal tensors associated with g^^. These are 
conveniently written in terms of a Lie-transported spa¬ 
tial triad (e^, , Cg) on 7 . The eigenvalues of the electric 

tidal tensor are 

Xf = (6) 

(with no summation over i); the magnetic tidal tensor 

has two equal but opposite eigenvalues, with magnitude 

(7) 

Here Rfiavp is the Riemann tensor of g^i,. 

From these scalars, we can obtain gauge-invariant per¬ 
turbative quantities by expanding g^y around g^y and 
7 around 70 . Of interest to us in this paper are the 
linear perturbations in ip, Af, and A®, which were de¬ 
rived in Refs. [52] and |5]. The explicit expressions in 
Schwarzschild coordinates are 


A 0 = 
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(9) 

( 10 ) 


( 11 ) 

























4 


AA® 


n\ro-M) »M( 5 ro - 9 M) ^ « 

(ro- 3 M) 2 (ro- 3 M )2 _ 3^)2 2 fl/\,{r, - iMf “ ^0 

[hfe,r - 2/ife,e - ^1^0,0] - + 1)2 3 [(^0 - ^M)h^^ ,. + 12(ro - 3 M)/i^^_J- 

0 2/0 Tq 2/g 


17 


Ofl/2 2/' QA/r\ 

2/0 r^iro-SM) 


ro ^t,eej 


2/„‘'V»"“«+2/V\3(.„-3M) 


(ro - M) 


'h4),ee^ 


( 12 ) 


where /o = 1 — ‘^M/tq, h^k = habk’^k’’, all functions are 
evaluated on the zeroth-order trajectory, and commas 
denote partial derivatives with respect to coordinates. 
The tidal invariants are related by AAf+AAf+AAf = 0. 


The quantities Atp, Af, and A®, as given above, are 
defined by performing an expansion at fixed 12. They are 
invariant under any smooth gauge transformation that 
leaves the helical Killing vector invariant; that is, for a 
smooth gauge generator the quantities are invariant 
if £^k^ = 0 (and hence = 0)- In practice, this 

means that to obtain the same result in any two gauges, 
we must ensure that the perturbation in each of these 
gauges possesses three properties: (i) helical symmetry, 
£khfj.v = 0, with k^ taking the form (ii) asymptotic 
flatness, (iii) regularity at 0 = 0 and tt, and (iv) con¬ 
tinuity of certain metric components. We explain these 
requirements in the next section. 


B. Computation of the retarded field 

We calculate the retarded field in a completed ra¬ 
diation gauge. After decomposing the perturbation into 
tensor harmonic i modes, we calcnlate the £ > 2 modes in 
an outgoing radiation gauge satisfying /iab^“ = 0, where 
2“ is the outgoing principal null vector. We calculate 
the £ = 0,1 modes, which physically correspond to the 
spacetime’s change in mass and angular momentum, in 
the asymptotically flat Zerilli gauge. We briefly review 
the reconstruction formalism here; more detailed descrip¬ 
tions are available in Refs. 0 [3SH371 ig. For the com¬ 
pletion, we make several points about continuity which 
have not previously been discussed in the literature. 

1. Reconstruction 

The £> 2 part of the metric perturbation is extracted 
from a Hertz potential, T, which in turn is calculated 
from the perturbed spin-2 Weyl scalar, tpo, by solving the 
separable Bardeen-Press-Teukolsky equation. The spin-2 
retarded Weyl scalar has the form 

'ipoix) = (13) 

with 


= dTTTriM* EE ^im [{£ - me +!){£ + 2)]i/2^H.,„(r<)Roo.,„ {r>hYm0, . (^a) 

0 e >2 m 

= STrimHu^Ao EE Atm - 1)(^+2)]i/22r,„(0, X 

i>2 m 

|[zmHro -f 2ro]i?H^,„(?'<)I?oof,„ (r>) + Ao[i?H^ _^(ro)i?oof,„(r')6»(r - tq) -f RHf.„(?')-R^,_„(?’o)6'(»’o - f)]}, (14b) 

EE 

i>2 m 

|[30rg - SOMrl -|-48M^rg - m^H^rg - 2Ao - 24Agr0(rg - M) + 6jmHrg(rg - M)]i?Hf,„ (?'<)Roof,„ (r>) 

-b2(6r^ - 20Mro -f 16M^r^ - 3roAo -f iml2Aoro)[RH^ _^(rg)i?oof.„(?’)0(T- - rg) -b (^o)-RHf.„ (?')6'(?'o - r)] 

+rlAl[R'^^mo)me,^{r) 9 {r - rg) -b i?((o^,„(»'o)I?Hf,„(T-) 6 »(rg - r) -b W[i?H^,„(r), i?oo,,„(r)] 5 (r - rg)] i. (14c) 
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Here A = - 2Mr\ Aq = Tq - 2Mro; sYt,m{d,(j)) 

are the spin-weighted spherical harmonics; the functions 
and Raoi^^{f) are the solutions of the homoge¬ 
nous radial Teukolsky equation that are ingoing at the 
future event horizon and outgoing at null infinity, and a 
prime denotes a derivative with respect to r; the Wron- 
skian, 


( 15 ) 


and the quantities given by 

~ A3W[i?H,,„(r-),i?oo,,„(r)]’ 

are constants, independent of r. 

The spin-weighted spherical harmonics sYt,miO, (j)) and 
their 9 derivatives at 0 = 7 r /2 are calculated analytically. 
On the other hand, we compute the functions „ and 
Root „ to more than 550 digits of accuracy using ex¬ 
pansions in terms of hypergeometric functions given in 
Ref. [48], namely 


OO 

^ anF{n + v+ 1 — ie,—n — v — ie,—l — 2ie]x), (17) 

n— — oo 

OO 

(-2z)^b^U{n + iy + 3-i€,2n + 2iy + 2;-2iz), (18) 

n——oo 


where a; = 1 — e = 2Mmft and z = —ex. We refer 
the reader to [451145] for the derivation of z/ (called the 
renormalized angular momentum), and the coefficients 
a„ and Here F and U are the hypergeometric and 
the (Tricomi’s) confluent hypergeometric functions. To 
minimize the computation time we hnd recurrence rela¬ 
tions between the different F's and C/’s, so that n=0 and 


1 are sufficient to find all other n’s in Eqs. (17 & 18). 


After obtaining the £m modes of i/'O) we calculate the 
radial part of the Hertz potential using the algebraic 


relation 


^ im — 8 


(-!)“(/ -k 2){e + !)£{£ - 1 )V'^,_„ -f 12 tTOMnV’fo 
[(/ + 2){e + l)£{i - 1)]2 -p 144m2M2fl2^ 

( 19 ) 

where the individual modes are dehned as the coeffi¬ 
cients in and 

V’o = J2e,m '4’tm{r) 2 yem{ 0 , (/))e"™^*. 

From the Hertz potential, the tetrad components of the 
reconstructed metric perturbation are obtained using 


r 2 

j^rec _ /•32,Tt I a2 


11 =y(g^vI; + g^Mr), 


urec _ 4 

^33 ~ ^ 


- 2fdtdr + 7292 


3{r-M)^ , /(3r-2M)„ , r'^-2M^ 
2r2 + 27^ + ^4 


h 


rec 
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dt - fdr 




( 20 ) 

( 21 ) 

( 22 ) 


where the subscript boldfaced numbers 1 and 3 respec¬ 
tively represent components along and m°‘ (the com¬ 
plex null vector on the 2-sphere). The other quantities 
appearing here are / = A/r^ and the operators 0 and 0 , 
the spin-raising and -lowering angular operators, which 
act on a spin-s quantity rj according to 


g?? = — {do + i CSC 9d^ — s cot 9) 77 , 

g?? = — (do — i CSC 9d^ + s cot 9) rj. (23) 


2. Completion 

As mentioned previously, the reconstructed metric per¬ 
turbation is incomplete. In the context of a Schwarzchild 
background, the reconstructed perturbation corresponds 
to the £ > 1 tensor-spherical-harmonic modes of the so¬ 
lution to the linearized Einstein equation with a point 
particle source. To complete the solution, we must add 
perturbations that satisfy the remaining, £ = Q and / = 1 , 
pieces of the Einstein equation. 

Appendix|^discusses a large range of gauge choices for 
the monopole solution. Here we add both monople and 
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dipole solutions in the asymptotically flat Zerilli gauge, 
as given in Eqs. (18)-(19) of Ref. [29]. The non-zero 
components of the £ = 0 mode are given by 


hSM _ 
'Ht < — 


2fiE 1 - 


2M 


hSM _ 

'ht > “ 

hSM Q 
^rr < 


, 5M ^ 
Urr > 


ro 1-^’ 
ro 

2iiE 


2^E 


and those of the £ = 1 mode are given by 


_ 2/rL^ sin^ 9 

— ys 

'o 

i^sj 2/rL^ sin^ 9 


(24) 


(25) 


The subscript < (>) indicates the metric perturbation 
in the region r < rg (r > rg). Outside the particle’s or¬ 
bit, at r > rg, the perturbation contains invariant shifts 
in the spacetime’s mass and angular momentum, corre¬ 
sponding to the contributions from the particle. Inside 
the particle’s orbit, at r < rg, the perturbations are pure 
gauge. We can see this explicitly by adding a gauge per¬ 


turbation £^ 5 ^,, with - r), 

which cancels the terms and given above. 

The reader may question why one choice of gauge is 
preferable to another, and why we cannot simply set the 
perturbation in the region r < rg to zero, given that we 
are calculating invariant quantities. The answer is that, 
as mentioned below Eq. (12), the quantities are invariant 
only within a certain class of gauges. As is well known, 
the quantities we are interested in (at least as they have 
been formulated historically) are invariant only within 
the class of gauges in which the helical symmetry is man¬ 
ifest, in the sense that the Killing vector takes the simple 
form ([^. It is also well known that these quantities are 
invariant only within gauges that agree on the time at 
infinity, which is usually resolved by working in gauges 
in which the asymptotic flatness is manifest, in the sense 
that limr_>.oo = 0(l/r) in asymptotically Cartesian 


background coordinates [50] However, here we point 
out that the quantities are invariant only within a class of 
gauges satisfying additional continuity conditions, which 
prove to be essential in our regularization scheme. 

Let us first restrict to the class of gauges with manifest 
helical symmetry. The gauge transformations that alter 


^ We stress “manifestness” here because this is only a matter of 
gauge choice; the physical solution will still possess a Killing 
symmetry even if we switch to a gauge in which that symmetry 
is not plain to see, and likewise for asymptotic flatness. 


our invariants even within this class are the ones which 
preserve the symmetry, in the sense that £khfj,u = 0, but 
which alter the Killing vector, in the sense that ^ 0. 
To our knowledge, the most general (continuous in r) 
gauge generator that accomplishes this is given by 

= {ait + a2(l>)t^ + {I3it +, (26) 

where ai and /3i are arbitrary constants. This generates 
a transformation with nonzero components 


Shtt = -2fai, (27) 

Sht^ = r'^Pisin^ 9 - fa2, (28) 

= 2r^P2 sin^ 9. (29) 

It alters the helical Killing vector to ^ + 6k^, 

where 

6ki^ = = -(ai -h na2)t>^ - (/3i -h (30) 


By rewriting k^ -|- Sk^ in the form ([^ (times an overall 
constant), we find that this shift simply corresponds to 
Q —>■ fl + 6fl + 0{fi^), with = (ai —/32 + f2ci2~/3i/f2)r2. 

One can reason that the transformations spoil invari¬ 
ance because our invariants are defined by an expansion 
at fixed 11, while these transformations alter il. This 
makes the transformations dangerous, because the gauge 
perturbation (27l~(29) is helically symmetric, satisfying 
£‘kShfj^i, = 0. So if two people compute one of the “invari¬ 
ants” in two different gauges related by this transforma¬ 
tion, they might each think they are working in suitable 
gauges. But the gauge perturbation contributes to the 
“invariant”, and the two people will obtain different val¬ 
ues from their calculations. To eliminate this possibility, 
we must impose conditions on that fix the freedom 
represented by the generator (26). First, we notice that 
the 02 term in Eq. (28) diverges (in regular Cartesian 
coordinates) at 0 = 0 and tt. So we can remove that 
freedom by requiring regularity at the poles. Similarly, 
the /32 term contains conical singularities at 0 = 0 and 
TT, with a deficit (or excess) angle of 27r,02; this is to be 
expected from the fact that in Eq. (26) it corresponds to 
the coordinate transformation (p ^ (j) — ^ 24 ’- 

This leaves us with the freedom in = aitt^ + fiitp^^. 
At this stage, we must consider jump discontinuities at 
r = rg. We replace a 1 with oi = ai+9'^ +ai-9~, and the 
analogue for /3i, with 9^ = 9[±{r — rg)]. This transfor¬ 
mation still preserves the helical symmetry of the metric 
perturbation, as the gauge perturbation it generates has 
the helically symmetric form 


Shft ——2f{ai+9^ + ai-9 ), (31) 

= r^(/3i+6»+-b/3i_6>“) sin^ 6». (32) 


As mentioned at the beginning of this discussion, we can 
partly abolish this freedom by imposing asymptotic flat¬ 
ness. Assuming we always work in an asymptotically flat 
gauge, the only remaining part of the unwanted freedom 
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is encoded in the gauge perturbations 


III. REGULARIZATION 


6htt =-2fai-9{ro - r), (33) 

Sht^ = r^/3i_6»(ro - r) sin^ . (34) 


How do we escape this last bit of freedom? Here we 
propose what we consider to be the most natural choice: 
we require continuity of the hu and ht^ components of 
the completion part of the metric. This brings us back 
to the £ = 0,1 solutions in Eqs. (24) and (25). We now 
see why we do not use the freedom in Eqs. (33|-(34) 
to eliminate and doing so would violate our 

continuity condition. 


In summary, the “invariant” quantities are invariant 
within the class of gauges satisfying conditions (i)-(iv) 
listed below Eq. (12), where we have now specified the 
continuity condition. We stress that this continuity con¬ 
dition (and the regularity conditions at the poles) are 
not imposed in a quest for regularity. Indeed, we hap¬ 
pily work with a perturbation that is not in a regular 
gauge: The rr component of Eq. (24) is plainly discon¬ 
tinuous across r = tq, and as we discuss in Appendix [X} 
this perturbation is also not regular at the horizon in 
horizon-penetrating background coordinates. The recon¬ 
structed perturbation is also discontinuous across r = tq, 
a fact which we find experimentally and discuss further 
in the next section. But none of these gauge irregularities 
affect the values of the invariants, because they are not 
generated by transformations of the form (26). 


Since the “invariants” do not have unique values, in 
principle we could equally well choose some other, less 
regular class of gauges to work in. Here we take a prag¬ 
matic stance: we want our invariants to agree with those 
that would be calculated by anyone working in a “nice” 
gauge that is asymptotically flat and continuous every¬ 
where away from the particle (at least for r > 2M, where 
the background coordinates are well defined). Clearly 
such a “nice” gauge falls within the class defined by our 
four conditions. 


Although we have specified the continuity condition 
at r = To, properly speaking it must be imposed at all 
radii. But it becomes important only at r = rg, because 
only there is one tempted to introduce a discontinuity. 
Especially in Kerr, setting the perturbation to zero in 
the region r < rg often appears the easy choice to make. 
It is doubly tempting because the perturbation is often 
already discontinuous across r = rg, as is the case for us 
here. 


If one is able to calculate the invariants from modes in 
the region r > rg, one need not worry about the issue of 
continuity. However, our regularization scheme will use 
the average of modes from r —>■ r()" and r r^. Hence, 
to ensure that our invariants have the same values as 
would be found in a “nice” gauge, we must impose the 
continuity conditions at r = rg. 


From the retarded field we now wish to calculate the 
invariants Atp, AAf, and AA®. In Sec. we defined 
these invariants in terms of the DW regular field, which 
was originally defined in the Lorenz gauge [S^. While 
it is easily generalized to any gauge smoothly related 
to Lorenz, it is less easily generalized to more singular 
gauges, among them the radiation gauge. In this section, 
we utilize the results of PMB to define a DW regular 
field in a completed radiation gauge; we show that the 
invariants constructed from this field are equal to those 
constructed from the regular field in any gauge smoothly 
related to Lorenz. Motivated by the parity structure of 
the transformation between the Lorenz gauge and radia¬ 
tion gauge, we compute the invariants using a mode-sum 
regularization formula that averages the £ modes from 
opposite sides of the particle. The validity of this proce¬ 
dure will be established in a companion paper. 


A. Detweiler-Whiting fields 

Radiation-gauge metric perturbations are known to 
have singularities extending away from the particle, ei¬ 
ther in the form of a string emanating from the particle 
or in the form of distributional singularities supported 
on a surface intersecting the particle dSl- In our calcu¬ 
lations in this paper, we work in the “no-string” gauge 
described by PMB; in the present context, this gauge is 
regular everywhere except at the sphere at r = rg, which 
supports a jump discontinuity and a delta function. 

Rather than working directly with the metric per¬ 
turbation, for our analysis we will instead refer to the 
gauge transformation from any gauge smoothly related 
to Lorenz. Let Xq be a point on the zeroth-order world¬ 
line 7 g. At a nearby point x^ = Xq + Sx^, the generator 
of the gauge transformation has the discontinuous form 

C'" = Kc+(a^o,^a;) -b Z!^(xo) + 0(elii£)]9+ 

+ [^o-(xo,^x) + Zt{xQ) -b O{e\n£)\0~. (35) 

Here e is a measure of distance from Xq, such that Sx^ ^ 
£. The quantities Z^(xo) are smooth functions of Xg. 
The other terms, are irregular at the particle. It is 
convenient to split them into pieces that become parallel 
and perpendicular to the worldline as e —> 0: 

^o± ='?||±(a^o,^a;)u'"(a:g)-bC±±(a;o,i5a:), (36) 

where U/^(xo)^±±(xo, 5x) = 0. The parallel piece behaves 
as 


C||±(a;o,(5a:) ~ ln(e). (37) 

Each of the perpendicular pieces and is bounded 

but has a direction-dependent limit 6x^ —)■ 0. The two 
are related by 

C±+(a;o,-fe) = -^±_ixo,Sx); 


(38) 















in other words, := has odd parity 

under the parity transformation that interchanges points 
on diametrically opposite sides of the particle. Its deriva¬ 
tives da^i3± have even parity; this is also true of 
even though has no definite parity. These general 
properties, derived by PMB, are all that we will re¬ 
quire. We refer the reader to Ref. [23] for explicit expres¬ 
sions. Given our construction of the metric, and the heli¬ 
cal symmetry, continuity, and asymptotic flatness of the 
(tweaked-for-asymptotic-flatness [50]) Lorenz-gauge solu¬ 
tion, it is clear that satisfies the condition = 0. 

The retarded field in the no-string gauge is given by 

hpiv = (39) 


choice of singular and regular fields in the new gauge: 

(41) 

= hl^. (42) 

With this choice, is guaranteed to satisfy properties 
(i) and (ii) mentioned above. 

Now let us consider our singular gauge transformation 
to the no-string gauge. Suppose we split into regular 
and singular pieces of the form 

=\{Zl + Z^) + 0{e), (43) 

es = e- (44) 


where we will use “L” to indicate that a quantity is in the 
Lorenz gauge. (An arbitrary gauge that is smoothly re¬ 
lated to Lorenz could equally well be used.) As promised, 
this gauge is discontinuous across the sphere at r = rg, 
due to the 9^ functions in and it contains a delta 
function on that sphere, due to derivatives of those func¬ 
tions. Our reconstruction and completion procedure can¬ 
not determine the correct S function on the sphere to 
make our computed metric perturbation a true solution 
to the linearized Einstein equation with a point-particle 
source, but PMB showed that the perturbation we obtain 
is pointwise correct away from r = rg. Hence, we shall 
always work away from the sphere. 

PMB established that in this irregular gauge, the linear 
correction Sz^ to the particle’s trajectory, relative to the 
trajectory in the Lorenz gauge, is given by the average 
of the regular part of the gauge transformation: 

=-i(Z()+ Z0). (40) 

This is the regular change in position induced by the 
singular generator 

How can we use this to define a Detweiler-Whiting reg¬ 
ular field in the no-string gauge? The DW regular field 
possesses three important properties: (i) the linear equa¬ 
tion of motion can be written as the geodesic equation 
in (ii) is a smooth vacuum solution to the 

linearized Einstein equation, and (iii) when evaluated at 
a point Xg on the worldline, and all its derivatives 
depend only on the causal past of Xg. To simplify the 
analysis, let us ignore the third property and focus on 
the first two. Consider a smooth transformation away 
from Lorenz. It corresponds to a coordinate transforma¬ 
tion x^ —>■ x^ — , altering the position of the particle 

to z>^ + 5z^^ with 5z^ = —^^{xg). It also alters the met¬ 
ric perturbation to -|- £^g^v. A priori, there is no 
unique way of splitting + £^gfiv into singular and 
regular pieces. However, if for a moment we ignore the 
singular field and consider a geodesic in a vacuum metric 
ff/ii/ + it is obvious that since the Einstein equation 
and geodesic equation are invariant, + £^gp,u will 
remain a vacuum perturbation, and z^ 5z^ will be a 
geodesic in gfj,ij + h^^ + £^gfj,,j. This motivates the natural 


where the 0(e) terms are smooth but for our purposes 
can be arbitrary. Then according to Eq. (40) we have 
Sz^ = —'Cij(a^o)- Hence, from the argument laid out in 
the case of a smooth transformation, the regular field 


( 45 ) 

will satisfy the desired properties (i) and (ii): the per¬ 
turbed motion is geodesic in g^i, -|- and is a 
smooth vacuum perturbation. With this choice, the sin¬ 
gular field is left to be 

(46) 


These will be our definitions of the singular and regular 
fields in the no-string gauge. 


B. Invariance of the invariants 


We now consider the construction of scalar quantities 
from (and its derivatives) on the worldline. Let A/ 
be one such quantity. It can be any one of Detweiler’s 
redshift, Aijj, AAf, AA®, etc., defined by substituting the 
no-string-gauge into Eqs. ([^~( 12). Before describing 


the computation of this quantity, we first confirm that it 
is equal to its value AI^ in the Lorenz gauge (and hence 
to its value in any gauge smoothly related to Lorenz). 

Let us begin by defining retarded and singular ver¬ 
sions of AI. We do this by extending AI away from the 
worldline; for example, in Eqs. ([^~( 12) we can multiply 
each term by an arbitrary function of t, r, 0, (j) that goes 
smoothly to 1 at a;o, and at points off Xg, we can let 
and its partial derivatives take their natural values. This 
defines a field AI[h^]. If we then replace the regular field 
with the retarded or singular field, we obtain fields AI[h] 
and A/[/i'^]. 

AI can be written as the limit of the difference between 
these two fields, 

AI{xg) = lim |A/[/i](a;) — A/[/i'^](a:)| . (47) 

X^Xq L J 


Note that while the extension is arbitrary, it must be the 
same in both terms on the right. Starting from Eq. (471, 
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we can establish our desired result, A/ = AI^, by trans¬ 
forming to the Lorenz gauge before taking the limit. In 
the extended fields A/, the only gauge-dependent quan¬ 
tities are the perturbations /i^^, and the only 
other dependencies are on the zeroth-order trajectory, 
which is gauge-invariant m, and on arbitrary smooth 
functions that dehne the extension, which we wish to 
keep the same in both gauges. Hence, the quantities 
constructed from no-string fields are related to quanti¬ 
ties constructed from Lorenz-gauge fields according to 


Ai[h] = Ai[h^] + Ai[£^,g^,] + A/[£5,5^,], (48) 

A/[h«] = Ai[h^^] + Ai[£^,g^,], (49) 


where we have appealed to Eqs. (391 and (46). Since 
AI is invariant under gauge transformations generated 
by smooth^ vectors satisfying £k^^ = 0, we have 
lima;_ya;Q A/[£j^ 5 ^i,] = 0. The terms in Eq. (47) can¬ 
cel, leaving us with 


A/(a;o) = lim \Ai[h^]{x) — Ai[h^^]{x)\ , (50) 

X^Xq L J 


which is our desired result: 


A/ = A/^. 


(51) 


We have now shown that the invariants constructed 
from the no-string-gauge regular field are equal to the 
ones constructed from the usual DW field in gauges 
smoothly related to Lorenz. However, we note that we 
could have split into and in any number of other 
ways and still obtained the same result. Eor example, we 
could have defined = 0, and the calculations in this 
section would have gone through just as well. In that 
case, the regular field in the no-string gauge would have 
been dehned to be identically equal to the Lorenz-gauge 
regular field. The message is that for the purposes of 
calculating gauge-invariant scalars on the worldline, in 
Eqs. ( |4^ and (46) the regular gauge vector is almost 
entirely arbitrary. So long as is smooth and helically 
symmetric, the choice has no effect on the invariants. 
Nevertheless, we think it is meaningful to define the reg¬ 
ular metric to be the one in which the perturbed motion 
is geodesic. 


C. Averaged mode-sum regularization 


Assured that we are calculating the correct quantities 
AI, we now describe our concrete method of performing 
that calculation. PMB showed that the gravitational self¬ 
force in the no-string gauge can be computed using mode- 
sum regularization in combination with averaging from 
opposite sides of the particle. Here we adopt the same 
method for our calculation of invariants. 


Eirst note that the direction of the limit in Eq. (47) 


is arbitrary so long as it avoids the sphere at r = rg. 


Therefore we can replace the limit with the average of 
two limits from opposite directions. 


AI = 


lim 
<5ai—>-0 


(A/[h]) - {Ai[h^^]) 


(52) 


where xq -I- Sx lies outside the sphere and xq — Sx 
lies inside it, and to keep the expressions compact, we 
have defined the average of a field F as {F){Sx) = 
5 [F(a:o -I- 6x) + F(xo — fe)]. 

Since we have a great deal of knowledge about the 
Lorenz-gauge singular field, we express in terms of 
it, giving us 


AI= lim 

<5x—>-0 - 


(AJ[h]) - (A/[h"^]) - (AJ[£5,ff]) . (53) 


Eollowing standard steps [46], we can rewrite this as a 
sum over scalar spherical harmonic modes. 


- OO 

= 9 E E \i^nh])i - {AI[h^^])i - {Ai[£^,g])i 

^ e=ok=± 

(54) 


where (A/)^(t, r, 0, ip) is the £th term in the spherical har¬ 
monic expansion of AI, summed over azimuthal number 
m, and (A/)^ = \im^^^±{AiY{t,6r,Tr/2,£lt) is its limit 
to the particle from r > rg or r < rg. (Because the orbit 
is circular, this quantity depends only on rg in the end.) 

For concreteness, let us examine Aip- We note that the 
singular held in the Lorenz gauge behaves as ^ 1/e. 
Since Aip^^ involves single derivatives of it behaves 

as Aiplh^^] ~ 1/e^. It is well known that this translates 
into the large-£ form 

{AY[h^Y)i = ^t±L + AtY + At\jL + 0 ( 1 /£ 2 ), 

(55) 


where L := £ -|- 1/2 and the regularization parameters 
are independent of I. From the general property 
^ Ine, we expect the same scalings for (A^/jT^gg])^, 
and we assume 


(A'!/’[£5sg])± — SAf^L + 5 Aq^ -|- 5AY_i±/L + 0(1/ L^). 

(56) 


Based on this, we write Eq. (54) as 
1 


A^=- 5] \(^Am)i-At,L-At-At,jL 


t,k=± 

1 

2 ' 


-l^Dt + Dt), 


(57) 


where and the D terms contain 

everything in the singular held modes that is not included 
in the A terms. Explicitly, 


Dt := Y.[{AY[h^^])i - Af^L - At^ - At.jL], 


e=o 


(58) 
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and the analog for 5d1^. Note that the parameters 
and depend on the choice of extension in defining the 
extended field and they must all be calculated 

in the same extension as is (A-0[/i])4- 

Using the method developed in Refs. [3SJE7j, we can 
numerically (and uniquely) determine the Ai± param¬ 
eters through large-£ fits of the modes of the retarded 
field. This allows us to numerically confirm the expected 
behavior (56). However, there is no obvious way of de¬ 


termining through numerical fits, and these D terms 
become the central concern. Here we eliminate that con¬ 


cern by contending that in the averaged form (57), the 


D terms cancel in the mode-sum formula, meaning they 
need not be determined. 

First we consider the contribution to d1^ from We 
argue that the parameters always vanish, regardless 
of the choice of extension. To see this, observe that 
has the schematic form | 3 ] 


h 


SL 


1 

P 


{dxy 

n3 


0{e), 


(59) 


where p ~ e is the leading-order geodesic distance from 
Xq to Xq + 6x^. Any quantity constructed from via 
single derivatives and multiplication by smooth functions 
will have the form 


SL 


s{x)dh^y 


6x 


{6x)^ {SxY 

p5 p7 


0(e), 


(60) 


where s{x) is smooth but arbitrary. Similarly, for two 
derivatives, as appear in the tidal invariants. 


s{x)d^h 


2,SL {5xf {5xf (fe) 


11 


+ - 


^ll 


-+0{e). (61) 


Now, since they are finite and hence can contain no posi¬ 
tive powers of L in their mode decomposition, the 0{£^) 
terms may seem likely to generate a D term. How¬ 
ever, these 0 (e*^) terms always have the discontinuous 
(direction-dependent) form (da:)"/p", with odd n. The 
sum of modes of a discontinuous quantity evaluated at 
the point of discontinuity, this case, con¬ 

verges to the average of the quantity on an infinitesimal 
circle around the point (51]. Because (da;)"/p” terms 
with odd n have odd parity around the point, they con¬ 
tribute nothing to the mode sum, and hence contribute 
no D. But the remaining terms in Eqs. (59)- (that do 
not vanish in the limit £ —>■ 0 ) only ever contribute non¬ 
negative powers of L to the mode sum [S21I13]; in fact, 
this should be true of any functions of the form {5x)'^/ff 
with odd p > 0 [S3]. So they too cannot produce a D. 

This argument should be easy to make precise and 
to extend to quantities constructed from any number of 
derivatives. Hence, we conclude that D terms vanish 
identically for any quantity constructed from the Lorenz- 
gauge singular field via the action of derivatives and mul¬ 
tiplication by smooth functions. This conclusion is inde¬ 
pendent of extension, since a change in extension merely 


corresponds to a change in the smooth functions multi¬ 
plying the derivatves. Therefore, we have = 0. 

Now let us turn to These D terms cancel in 

Eq. (57) if they have the symmetry property SD't = 
—5D'^. Establishing this property requires a detailed lo¬ 
cal analysis of A'i/;[£jgp], which will be carried out in a 
companion paper. In the meantime, we present a plau¬ 
sibility argume nt. Fro m the general properties of de¬ 
scribed in Sec. 
form 


HI A 


it follows that A^/>[£jgp] has the 




' [1 -I- i5x -I- ((5x)^]9^^s ■ 
+ Cs + o(£°) 


[l-b(5a;]9^s 


(62) 


where “o(e°)” means “goes to zero in the limit e —>■ 
0”. Now let us assume that all positive powers of L 
in (A' 0 [£jgp])^ arise from the negative powers of £ in 


Eq. 


and that SOf arises from the order-£° terms. 


Eor the purpose of this sketch, consider only the most sin¬ 
gular terms in d^g. Given this restriction, order-£*^ terms 
in Eq. (62) can only arise in the forms {6x)^d'^^s, Sxd^s, 
and ^s- Recall that the most singular terms in d^s have 
even parity, meaning 6xd^s has odd parity; hence, this 
term will vanish upon averaging. Similarly, the most sin¬ 
gular term in has odd parity, and so {6x)^d'^^s also 
vanishes upon averaging. Last, the ^5 term in Eq. (62) 
can only come in the form ^gd^gap, stemming from such 
a term in Since u^d^ga /3 = 0, only con¬ 

tributes to this term; and since has odd parity, this 
contribution also vanishes under averaging. Therefore, 
if 6D^ is generated solely from the terms in Eq. (62) 
that have finite, nonzero limits as £ —>■ 0 , then averaging 


eliminates the D terms in Eq. (57), as desired. 


Following these arguments, we now allow ourselves to 
set D± to zero in Eq. (57). This gives us our final formula 
for A-tp: 


lip 


^ e,k=± 


(63) 


Here we have used = 0, a fact which we find nu¬ 

merically. 

Following identical arguments for the tidal invariants, 
we obtain 


^-^—2 ~ A 2 kL ^ — ~ ^Ok 

^ i,k=± ^ 


, (64) 


where A A is any of AAf or AA®. 

When implementing these formulas, we first determine 
the A parameters numerically. We find that the param¬ 
eters, with our choice of extension, agree precisely with 
the analytically derived parameters in the Lorenz gauge 
in Ref. | 8 ] , with the lone exception of A^^. After making 
this determination, in practice we use the analytical pa¬ 
rameters except in the case of which we determine 
numerically. 
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As an aside, we note that one need not always define 
the modes in these sums as scalar spherical harmonic 
modes. It is computationally simpler to leave the har¬ 
monics in the “natural” form they arrive in from the re¬ 
construction procedure, which is a mixture of harmonics 
of various spin weights; using this “natural” form avoids 
having to re-expand spin-weighted harmonics into scalar 
harmonics. Although strictly speaking the arguments 
made in this section assume the modes are defined in 
the scalar harmonic basis, we have experimentally found 
that our results are unaltered by changing between the 
“natural” and scalar basis. It is has been suggested to us 
that this may hinge on our averaging from two sides of 
the particle [53j , as such averaging is known to annihilate 
D terms that arise in some cases when nonscalar bases 
are used. However, we leave exploration of this issue to 
future work. 


IV. EXTRACTION OF PN COEFFICIENTS 


In our concrete implementation of the mode-sum for¬ 
mulas (63|-(64), we follow the method of Ref. [35] ■ We 
fit the i modes to polynomials in L to numerically de¬ 
termine the regularization parameters Ai±, and to gain 
faster convergence, we also fit the power-law “tail” of the 
singular field, which is made up of negative powers of L. 

Because we wish to extract pN coefficients, we calcu¬ 
late the invariants at very large orbital radii. Specifically, 
we calculate the invariants for 144 different values of H 
by placing the particle at 


R = 1, 2, 3, ..., 9 X 101849,20.....33^ (gg) 

where R = is ^ coordinate-independent mea¬ 

sure of orbital radius, which in Schwarzschild coordinates 
is equal to the orbital radius rg of our zeroth-order trajec¬ 
tory. To obtain our desired accuracies of 1 part in 10^°°, 
we calculate £ modes up to a maximum ^rnax = 200 . 
These extreme accuracies are made possible with Math- 
ematica’s arbitrary-precision algebra. 

To extract the pN coefficients from these results, we fit 
the high-precision numerical data to the pN series in the 
following manner: First, we utilze the known pN series 
calculated by Bini and Damour, which were obtained to 
8.5pN order for Sij} and 7.5pN order for the tidal invari¬ 
ants. We subtract these series from the numerical data 
and then fit the result to a pN series j/^ log*(j/), where 
y = 1/i?; i/2 is the pN order of a term. In perform¬ 
ing the fit we first extract as many analytical coefficients 
of the highest power of log(j/) for a given pN order as 
possible. This is done before anything else because the 
analytical coefficient multiplying the highest power of log 
is always either a rational number (for integer powers of y 
in Stjj and AA^, and for half-integer powers of y in AA^) 
or a rational number times tt (for half-integer powers of 
y in 6't/’ and AA^, and for integer powers of y in AA^). 
Extracting two whole numbers (p and q oi p/q or pir/q) 


from the numerical value requires less precision than ex¬ 
tracting rational coefficients of transcendentals like tt^, 7 , 
log(2), log(3), etc., which make up the pN coefficients of 
less than maximum power of log( 2 /) for a given pN order. 
Once the coefficients of the highest power of log(p) are 
calculated we subtract them (along with the known pN 
series) from the numerical data and fit the resultant to 
the pN series (without the terms whose coefficients are 
known) and extract the coefficients of the second high¬ 
est power of log(y). These coefficients usually have tt^, 
7 , log(2) and log(3) in them. We repeat these two steps 
until as many analytical coefficients are extracted as pos¬ 
sible. The third step is to subtract all the analytical coef¬ 
ficients from the data and fit the remaining pN series, the 
coefficients in which are determined as numbers with fi¬ 
nite accuracy rather than as analytical expressions. This 
procedure, explained in more detail in Ref. |20j . enables 
us to obtain the pN expansion to significantly higher or¬ 
der than if we performed a straightforward numerical fit 
of the data to the series y^ log*(j/). 

In Figs. 0-§- we compare the pN-approximated in¬ 
variants with their exact (up to numerical error) values 
computed in Ref. [ 8 ] , and we find good agreement all the 
way to near R = 1/y = AM. More precisely, we find 
that our pN expressions are accurate to more than 4 dig¬ 
its up to the innermost stable circular orbit (i? = 6 M), 
and the approximation then worsens to less than 1 digit 
at i? = AM. We also examine the efficacy of exponen¬ 
tially resumming the pN series (as described in Ref. |54 ) ). 
Figs. (5]|9) show the relative errors in the pN series and 
in the exponentially resummed series for the tidal invari¬ 
ants at i? = 4, 5, 6M. We see that resummation does not 
markedly improve the accuracy of the series. As a consis¬ 
tency check, we have verified that the sum of the electric- 
type tidal invariants, AAf, AXf and AAf', is zero for any 
pN order. We also find that the contributions to this sum 
from the retarded field, prior to regularization, vanishes 
mode-by-mode for each £. 

Because of the great length of the final expressions, we 
place the analytical coefficients in Appendix and the 
numerical coefficients in the supplementary text files. 


V. DISCUSSION 


In this article, we have developed and implemented 
tools for computing invariant quantities in a “no-string” 
completed radiation gauge. At the level of formalism, 
we have defined a Detweiler-Whiting regular field in the 
no-string gauge, and we have shown that the invariants 
constructed from this field take the same value as they 
would if constructed from the DW field in any gauge 
smoothly related to Lorenz. At the level of practical 
calculation, we have described a simple mode-sum regu¬ 
larization scheme involving averaging from two sides of 
the particle. In the cases of the spin-precession and tidal 
invariants, the scheme is encapsulated in Eqs. 


The validity of this scheme will be rigorously established 
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in a companion paper. We have also identified continuity 
conditions that must be imposed on the completion part 
of the metric perturbation in order to obtain the correct 
values for the invariants, and we expect that an extension 
of these conditions will be vital to calculations in a Kerr 
background. 

By applying our regularization scheme, we obtained 
our main results: very high order pN expansions of 
the linear-in-mass-ratio corrections to the spin-precession 
and tidal invariants for a particle in circular orbit about 
a Schwarzschild black hole. These expansions were ex¬ 
tracted from numerical results at very large orbital radii, 
from 10^®M to 9 x Using the analytical solution 

to the Teukolsky equation found in Ref. [48] , we were able 
to obtain these numerical results with accuracy greater 
than 1 part in 10^°°. This data was then used to nu¬ 
merically fit and extract the analytical and numerical 
coefficients in the pN expansion of the invariants. 

To increase the accuracy of these numerically calcu¬ 
lated quantities (at such large radii), one would have to 
drastically increase the maximum computed i mode.This 
difficulty arises because of the large-£ behavior of the sin¬ 
gular field, which makes the mode sum converge slowly 
unless one subtracts both the positive-power-of-L terms 
and a large number of tail terms. The situation here con¬ 
trasts with calculations of fluxes, in which the singular 
field plays no role. To compute fluxes with an accuracy 
of 1 part in 10®°° in Ref. [T9| , going to £max = 40 was suf¬ 
ficient, whereas to compute conservative invariants here, 
we went to £max = 200. 

Instead of hunting greater accuracy, the next step will 
be to calculate the conservative invariants, along with the 
pN expansion of the Detweiler redshift, in Kerr space- 
time. This should provide a substantial improvement 
to current EOB models. A computation of spin preces¬ 
sion for a particle in circular orbit about a Kerr black 
hole will provide the linear-in-mass-ratio, spin-spin part 
of the effective Hamiltonian and improve the EOB model 
for tidally interacting spinning binary systems. 
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Appendix A: Metric perturbation corresponding to 
the change in mass of the spacetime 


The 1 = ^ perturbation we used in Sec. |IIB^ was suit¬ 
able for our purposes, but by transforming to ingoing 
Eddington-Finkelstein (EF) coordinates, one can see that 
it is irregular at the future horizon. In addition, it is dis¬ 
continuous across r = tq. Such gauge singularities do not 
affect the calculation of invariants, but it is worthwhile 
to construct the perturbation in a gauge that is mani¬ 
festly regular. Since, it is easy to check the regularity of 
the metric at the future event horizon in the ingoing EF 
coordinates {v,r,9,4>), we start with the Schwarzschild 
metric in the those coordinates: 


ds^ = — _|_ 2dvdr + r'^dO^ + sin^ 9d(j>^. 

(Al) 


We know from the solution in Sec. dIB 21 that across 
r = ro, there is an invariant shift of the spacetime’s mass 
by an amount fj.E, and everything other than the mass in 
the perturbation is gauge. (The jump in mass can also 
be discovered from the invariant definitions of mass in 
Ref. |55j.l So we begin with a straightforward variation 
of the mass in Eq. (Al), 


,ef at \ 

Kv = -6»(r-ro), 


(A2) 


with the rest of the components being zero. Note that 
the linearized Einstein equation only teds us about the 
jump in mass across the particle’s orbit; it does not ted 
us anything about the total mass, which can be freely 
altered by adding a trivial perturbation of M everywhere 
in spacetime. Here we have chosen to set the shift in mass 
to be zero at r < tq. 

We now perform a gauge transformation, —>■ hab = 

so that hab is continuous across r = tq, 
regular at the horizon, and asymptotically fiat at spatial 
infinity. Our strategy is similar to that of Ref. |3D]. We 
choose a spherically symmetric gauge vector (which re¬ 
spects the symmetry of the problem), 


C = {P{v, r), Q{v, r), 0,0). (A3) 


The term, Eab := + '^bia, is 
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^ (2g- ^g-2/p) 

(P + Q'-/P') 

0 

V 0 


[P + Q'- fP') 
2P' 

0 

0 


0 0 ^ 

0 0 

2rQ 0 

0 2rQ sin^ 9 j 


where a dot and a prime denote derivatives with respect the form of Eq we set Q = Q{r) and P = P{v). 

to V and r, respectively, and / = (1 —2M/r). To preserve Using this, we now have. 


-•ab — 


^(-^Q-2/p) (p + Q') 0 

{P + Q' 


To make the metric perturbation independent of the v 
coordinate, we have P = av where a is a constant which 
needs to be determined from the three conditions men¬ 
tioned earlier. We now split our gauge vector in half, 
5 = ^<0(ro — r) -I- ^>9{r — tq). Since we want a regular 
perturbation at the horizon, we choose Q on the inside 
to have the form (r — 2M)*, where i is a positive integer 
which will be fixed later, and to guarantee flatness at spa¬ 
tial infinity we choose Q to have the form l/(r — 2My, 
with j a positive integer which will be fixed later. With¬ 
out loss of generality, we write 

C<= (a<v,l3^( , 0 ,oV 


= oi>v,y> 


tq — 2M J 

ro - 2M^^ 
r-2M 


,0,0 , (A4) 


in which and need to be calculated. Imposing 
continuity of hgg or at r = tq tells us that 


/3< = /?> = /3. 


(A5) 


Imposing the asymptotic flatness condition gives us a> = 
0. Imposing continuity of hyy fixes a<. 


a< = 


—fiE 


rofo ’ 

and imposing continuity of hyr gives us 

liE 


/3= T 


1 + J 


(A6) 


(A7) 


We have thus calculated the metric perturbation which 
is regular at the horizon, continuous across r = ro and 


0 

2rQ 


0 2rQ sin 9 j 


flat at spatial infinity. We now write this metric per¬ 
turbation in Schwarzschild coordinates {t,r,9^(j)). For 
r < ro. 


htt — 2P/i 
htr = P/r 


1 _ 2M 
r 

M 

/ r - 2M y 

To — 2M 

j) 

\ro-2M) 


r-2M 


{ir-2M) 

ro — 2M r{i + j)(r — 2M) \ro — 2M^ 

2En{ir-M) fr-2MV 
{i+j){r-2MY \ro-2M) ’ 

2rE^x / r - 2M y 
(i+j)\ro-2M^ 




2rEii sin^ 9 / r — 2M \ * 
(*+j) Vro-2M, 


For this metric perturbation to be regular 
event horizon, i > 2. For r > ro, 


liO llictxiv.., pCJ-tUiUCXl 

horizon, i > 2. For r > ro, 
htt = 2.Eii 

/ / ' 

(j r + 2M) 


(AS) 

at the future 


htr = Efl 

hrr = 2EfI 


I M /'ro-2My 

r r^(i+j) \r — 2M J 

2 


\o - 2My 


r — 2M r{i+j){r — 2M) \r — 2M 

U r + M) _ 

2M)2 \ r-2M 


, yj , -r I ^0 - 2M\^ 

(r — 2M)2 (*+j)(?' — 


2rEfi /ro-2My 
(* + j) {r-2M J 

2rP/i sin^ 9 f tq — 2M \ 
(* + j) \r-2M J 


(A9) 
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For this metric perturbation, j > 2. We find that this 
metric perturbation is not manifestly static, having a 
nonzero time-independent htr component. Of course, 
this is simply a gauge artefact, as we know from the so¬ 
lution in Sec. |IIB 2 that the metric is static; in the gauge 
of this section, the timelike, hypersurface-orthogonal 
Killing vector is simply no longer but some appro¬ 
priate -I- St^^. 

The perturbation we have found is actually a large 
class of (physically equivalent) perturbations, corre¬ 
sponding to a variety of gauges labelled by i and j. How¬ 
ever, no two members of this class are related by the he¬ 


lically asymmetric gauge vectors discussed in Sec. IIB2 


and in line with that discussion, we have found exper¬ 
imentally that the final result for any of the invariants 
is independent of the values of i and j. Furthermore, 
again in line with Sec. IIB 2 we have also verified ex¬ 


perimentally that continuity of the components of the 
£ = 0 solution need not be imposed except on (or 
htt)', discontinuity of other components does not affect 
the invariants. 

It is not difficult to see that using the spherically sym¬ 
metric gauge vector given in Eq. (A3), it is not possible 
to find a gauge in which the metric perturbation is man¬ 
ifestly continuous, regular at the horizon, asymptotically 
flat, and static. This agrees with the extensive analy¬ 
sis (focused on the Lorenz gauge) by Dolan and Barack 
in Ref. [5S]. Perhaps surprisingly, this appears to stem 


from our choice of total mass of the spacetime. There 
do exist £ = 0 solutions with all the above desired prop¬ 
erties, but with a total mass of M -|- dM -|- p,E rather 
than M -|- p.E. These slututions, first discovered in the 
Lorenz gauge by Berndtson [SS], correspond to having 
transplanted a small part of the black hole’s mass into 
the perturbation, such that the mass of the background 
spacetime, M, differs from the black hole’s physical mass, 
Mbh = M + dM. While one can use these solutions, one 
must be mindful that Mbh ^ M in order to avoid any 
hiccups. 


Appendix B: Analytical results 


Most terms in our pN expansions are obtained as real 
numbers with a finite numerical precision. However, 
many terms we find in exact, analytical form. In this 
appendix we present those analytical coefficients. The 
numerical coefficients are presented in the accompanying 
text files. 

We present the analytical coefficients of 6ijj for the fol¬ 
lowing pN orders: (i) all the coefficients that were already 
known, which go up to [2^ (as later corrected in 

Ref. [30]), (ii) the highest power of log( 2 /) for each pN or¬ 
der from J/® to and (iii) the second-highest power 

of log(y) for each integer-pN-order from y® to y^^. 


Aijj = y'^ - 3y^ - y y'^ 


/-6277 4961og(2) 2047l7r2 , A . f 87055 6536297r2 527 


37721og(2) 7291og(3) 261og(y)\ g 265367 r ^ 3/2 / 149628163 76287 2977619477r2 14079877r4 

105 14 5 1575 ^ V 18900 21 393216 524288 

45561og(2) 12879 log(3) 38141og(y)\ 7 11341l7r75/2 /403109158099 749094627 34247® 

21 35 21 ~ 22050 ^ V 9922500 70875 25 


1646739794577r® 1609347643177r4 340681718 log(2) 8696967 log(2) 58208 log^ (2) 199989 log(3) 

353894400 335544320 1819125 1575 105 352 
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9765625 log(5) 1344C(3) 37454731 log(y) 434848 log(2) log( 2 /) 3424 


28512 

11795912067r 17/3 

3274425 ^ 

844, 2/ q 

-^log {y)j + 


5 70875 

/ 2227389947 log(j/) 
~ 606375 

232643687 r log(?/) ;^g /2 
■ 165375 


3376 


7 log(y) 


1575 25 

1538864 log(2)log(y) 


7log(2/) + 


856 log'(y) 


25 


y 


28431 


15 ' 11025 49 

5724 079403437 l og(y) , 355702967 log(y) 
496621125 


log(3) log(y) 

239211992 log(2)log(y) 


54832464 log(3) log(y) 8892574 log' (y) 


13475 

26744864 log(2) log'(y) 


11025 
3114128 log^(y) 


23625 

6900152727 log'(y) 
385875 


33075 

14320872 log(2) log'(y) 
42875 


10 231620397r log(y) 31/2 

^ 128625 ^ 


11 72162982310387 rlog(y) 23/2 

^ 2269176525 ^ 


11025 654885 

^ 1059800138707 log' (y) , 


62282567 log'(y) 


756392175 11025 

/30579194748876386 log' (y) 
1491353238375 


1108809 2 , , 345007636log3(y) 

log(3) log (y) - 


343 


1157625 


,12 


y 


14308502247r log'(y) 25/2 ( 58168822370617659002 log'(y) 489023054528727 log'(y) 

2480625 ^ V 2214659558986875 3781960875 


18295600819864 log(2) log' (y) 
1620840375 

18776756797975 log'(y) 27/2 
121550625 ^ 


116958256614 log(3) lo g'(y) 24451152726436 log^(y) 


5187875 
1857799285622072434log^ (y) 
33219893384803125 

,4 


,13 


11345882625 
267353384327 log^(y) 
17364375 


160490433344log(2)log^(y) 3341917304log^(y) 


52093125 17364375 

100052121353936113534941og^(y) 1610456327 log^(y) 

+ 1 , + 19845 


14 , 50447196163003061277 log'(y) 29/2 

^ ■ 39313483295625 ^ 


143952871334146875 
28829034 log(3)log^(y) 20130704 log'‘(y) 


1537620535168 log(2) log^(y) 
364651875 


2401 

27617258524023877 log^(y) 33/2 
38288446875 ^ 


^ 15 858623818390477 log^(y) ^ 31/2 131271672785838026 log'‘(y) ^ ig 

^ 5469778125 ^ 39313483295625 ^ 


19845 

27264557008 log®(y) 17 


86821875 


-y 


225632261039219879361277 log®(y) 35/2 
' 81732731771604375 


23312186833892 log®(y) ^ ig 183812819618619277 log^(y) 37/2 
9116296875 ^ 574326703125 ^ 


(Bl) 


For AA®, we present (i) all the pN coefficients (known 
from (31]) up to y'®/', (ii) the highest power of log(y) 


for each pN order from y^® to y'*®/', and (iii) the second- 


highest power of log(y) for each pN order from y^ 
yi7. 


to 


AAf = 2s= + 2s‘+ 


71779 71977^ 15361og(2) 7687 384 log(y) ^7 


227 59377 " 

^3 ^ y " 'V 4800 256 ' 5 ' 5 ' 5 J 


54356247 

2835 


^5629703 1008787772 52481og(2) 171527 29161og(3) 85761og(y)\ g 273927747/2 

100800 24576 7 105 7 105 ^ 175 ^ 

4692901483772 8774321og(2) 20898log(3) , 219337377^ 6746904013 2717812 log(y) ^,9 25411677 39/2 

7077888 2835 

/582414031287 876544 

' 5457375 17^ 


y 


-y 


7 log( 2 ) 


7 1048576 7257600 2835 J " 1225 

113134518813241772 6396680456 log(2) 6028101 log(3) 665335740577^ 


12288C(3) 97656251og(5) 876544 21913672 


5 

219136 

175 


-f 


3168 


175 


19818086400 
:iog'(2)- 


5457375 1120 67108864 

1964481413350639 29120701564 log(y) 


, , , 438272, , 

7 log(y)-log(2) log(y) - 


175 48898080000 5457375 

54784 2 / aV.io 597703934677 21/2 /170466930577log(y) 

{.y) j y ^ ^ 1 oAcri/icA 


175 


3274425 




9459450 
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39488 


7 log(j/) 


18028096 log(2) log(j/) 227448 


315 
4689510477r 


log(3)log(y) 


98721og2(y) 


+ 


3675 49 

4382727r3 937902087r log(2) 234475527r 


18375 525 

1699541000487 log(j/) 
9823275 

/4255324947293257197r 
' 14913532383750 

2139283767rlog(?/)\ 25/2 
385875 J ^ 

46895104 log^(y) 

55125 

8024633698757671 log(2) 
3094331625 
10566741767 log^(?/) 


18375 

30332610032 log(2)log(?/) 


1964655 

42785675277r 


18375 

1630044 

49 


log(y) \ 




315 

,23/2 


y 


11 


/ 642655659828430971 


+ - 


524431908000 
53463358651678127 log(j/) 


log(3) log(j/) + 


884978844750 
42488525012 log^(i/) 


13 


385875 11025 

815629077100102 log2(?/) 
49165491375 

42557714670757025410503171 

3212732800236960000 


16089104713 1256184276871 log(2) 


9823275 

2956824 


1157625 

937902087log^(^) 

’ 18375 ' 

928193260264771 


log(3) 


343 

187580416 log(2)log^(?/) 
18375 

308341847i3 


28014525 


21190572 


log(3) 


46409663013271 log(y) ^^ 27/2 


343 

1684425088 log(2)log^(7/) 


28014525 
8870472 


2695 

235234718516611549 log2(?/) 


165375 128625 343 

484644167389753134815371 log(?/) 40142209024771 log(y) 


log(3)log {y) 


2130504626250 
528337088 log^ ( 7 /)' 


.14 


y 


496125 

375160832713 log(y) 8028441804871 log(2) \og{y) 


59057588239650000 
2 


1929375 


55125 


1929375 


1003555225671 log^(?/)\ 29/2 /255969177552685551119891og^(7/) 3268309988247527 log^(j/) 


1929375 

5873797739780528 log(2) log^(?/) 


387249042885705000 


4862521125 


63571716 


log(3) \og^{y) - 


163415499412376 log3(j/) 


+ - 


34037647875 343 

4192025770928737878384771 log(?/) 674951603104771 \og{y) 


,15 


14587563375 
9936857824713 


111963344371003125 40516875 385875 

21798132745671 log(2) log(y) 23063227271 log(3) log(j/) 16873790077671 \og^{y) 


3472875 


2401 


40516875 


,31/2 


/ 681852850098374173148 log3(y) 802844180487 log3(y) 160568836096 log(2) log3(?/) 

33219893384803125 5788125 ' 5788m 

10035552256 log^ ( 2 /)' 




2/36 + 


/ 5777417607325841082973552371 log( 2 /) 1472383886624419184771 log{y) 


5788125 

247867707836464713 log ( 2 /) 

4862521125 ' 

36809597165610479671 log^ {y) 
5616211899375 
1195410487168 log(2) log3(2/) 


V 124777599295357008000 
7139241235001082414471 log(2) log( 2 /) 


5616211899375 
165286461671 log(3) log( 2 /) 


2/33/2 ^ 


117940449886875 

^ 1170236981712855093944183 log6(2/) 


3023010298017084375 
230632272 log(3) log^ (y) 393054276976 log'‘ (y) 


2401 

31444342158087 log6 ( 2 /) 
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121550625 2401 72930375 

00^^01 8590432731136 log® (/z) 


17 859043273113671 log^( 2 /) 35/2 

^ ~ 607753125 ^ 


4868681792431448132log^( 2 /) ig 4997543267316871 log^( 2 /) 37/2 
353821349660625 ^ 1418090625 ^ 


3038765625 


-y 
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(B2) 


For AAf, we present (i) all the pN coefficients (known 
from |3T]) up to 2 /^®^^) (ii) the highest power of \og{y) 


for each pN order from y^^ to and (iii) the second- 

highest power of log( 2 /) for each pN order from to 
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AAf = -y^ - - —y^ 
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For AAf, we present (i) all the pN coefficients (known 
from [31]) up to (ii) the highest power of log(y) 


for each pN order from to and (iii) the second- 

highest power of log(y) for each pN order from y^° to 
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FIG. 1: Fractional difference between the different pN-order-approximations and the numerical value of Sip. 
Color-code is as follows - Black: 19pN, Magenta: 18pN, Blue: 17pN, Green: 16pN , Orange: 15pN, Red: 14pN. 
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FIG. 2: Fractional difference between the different pN-order-approximations and the numerical values of AAf, AAf", 
and AAf. Color-code is as follows - Black: 20pN, Magenta: 19pN, Blue: 18pN, Green: 17pN , Orange: 16pN, Red: 

15pN. 
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FIG. 3: Fractional difference between the different pN-order-approximations and the numerical value of AA^. 
Color-code is as follows - Black: 20.5pN, Magenta: 19.5pN, Blue: 18.5pN, Green: 17.5pN , Orange: 16.5pN, Red: 

15.5pN. 








FIG. 4: Fractional difference between the numerically calculated value and Taylor-summed pN series (represented by 
red line), and the fractional difference between the numerically calculated value and exponentially-resummed pN 
series (represented by blue line), for different pN-order. These plots are at i? = 6M. 
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FIG. 5: Fractional difference between the numerically calculated value and Taylor-summed pN series (represented by 
red line), and the fractional difference between the numerically calculated value and exponentially-resummed pN 
series (represented by blue line), for different pN-order. These plots are at i? = 6M. 









FIG. 6: Fractional difference between the numerically calculated value and Taylor-summed pN series (represented by 
red line), and the fractional difference between the numerically calculated value and exponentially-resummed pN 
series (represented by blue line), for different pN-order. These plots are at i? = 5M. 
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FIG. 7: Fractional difference between the numerically calculated value and Taylor-summed pN series (represented by 
red line), and the fractional difference between the numerically calculated value and exponentially-resummed pN 
series (represented by blue line), for different pN-order. These plots are at i? = 5M. 
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FIG. 8: Fractional difference between the numerically calculated value and Taylor-summed pN series (represented by 
red line), and the fractional difference between the numerically calculated value and exponentially-resummed pN 
series (represented by blue line), for different pN-order. These plots are at i? = AM. 
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FIG. 9: Fractional difference between the numerically calculated value and Taylor-summed pN series (represented by 
red line), and the fractional difference between the numerically calculated value and exponentially-resummed pN 
series (represented by blue line), for different pN-order. These plots are at i? = AM. 



